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PREFACE 


This  report  was  prepared  by  Dr.  Theodore  S.  Bolls,  of  the  State 
University  of  New  York  at  Oneonta,  In  residence  at  the  Rome  Air  Development 
Center  (RADC),  under  the  1978  USAF-ASEE  Summer  Facility  Research  Program. 
The  work  presented  Is  a part  of  an  RAOC  program  to  develop  Bayesian  and 
other  statistical  techniques  based  on  the  use  of  prior  data  for  practical 
application  to  Reliability  Demonstration. 


1.  INTRODUCTION 

1.1  We  consider  equipment  with  exponentially  distributed  tlme-to-fallure, 
l.e.,  the  probability  density  function  of  the  tlme-to-fallure  Is  given  by 

(1.1.1)  ( 1 1 0 ) - P'1  exp  (- 1/6 ) . t s o;  0 s o, 

where  the  parameter  0 Is  the  mean-time- to- fai lure  (MTTF)  of  the  equipment. 

We  assume  that  0 itself  has  a prior  distribution  of  the  inverted  gamma  type, 
i.e.,  the  prior  probability  density  function  of  0 Is  given  by 

(1.1.2)  g (0;  \,y)  * -— ) exp  (-y/P);  P •*  o;  X o,  y > o, 

where  \ is  the  shape  pa rameter  and  > is  the  scale  parameter  of  the  prior 
distribution. 

1.2  Bayesian  Reliability  Test  Plans  based  on  the  prior  (1.1.2)  have  been 
developed  by  Schafer  et  al  [4]  and  Goel  [1]  under  various  combinations  of 
risks.  The  implementation  of  these  plans  require  the  estimated  values  of 

\ and  y in  (1.1.2).  Since  the  true  MTTF  P of  an  equipment  is  not  observable, 
we  cannot  directly  fit  existing  data  to  the  inverted  gamma  distribution 

(1.1.2) .  To  get  around  this  difficulty,  we  consider  the  probability 
function  of  the  number  of  failures  r in  a fixed  T,  given  0.  Because 
of  the  exponential i ty  assumption  (1.1.1),  this  probability  function  is 
Poisson  with  parameter  T/0,  1.  e. 

(1 .2.1)  PT  (r(p)  - -J—  (T/p)r  exp  (-T/0),  r * o,  1 . 2,  . . .;  T " o. 

T r j 

Thus,  the  unconditional  probability  function  of  the  number  of  failures  r 
In  a fixed  time  T Is 

(1.2.2)  PT  (r)  » / Pt  (r|0)  g (0;  \,  y)  dP 

1 o 
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By  using  (1.1.2)  and  (1.2.1)  and  performing  the  integration,  we  obtain 
(1.2.3)  PT  (r)  * ( X r ) ( — )X  (— — — )r  , r * 0,  1,  2,  which 

is  a negative  binomial  distribution  with  parameters  \ and  T/(T+>).  If 
existing  data  on  a type  of  equipment  are  of  the  form  "number  of  failures  in 
a fixed  common  time  T",  then  the  parameters  \ and  y can  be  estimated 
by  using  (1.2.3).  Schafer  et  al  [3]  used  the  method  of  moments  for  this 
purpose,  whereas  Goel  and  Joglekar  [2]  used  the  maximum  likelihood  method. 
1.3  An  extreme  and  rather  hypothetical  case  results  when  we  keep  the 
number  of  failures  r fixed  and  observe  the  time  T until  the  rth  failure. 
Since  T is  the  sum  of  r exponential  random  variables,  its  probability 
density  is  gamma  with  parameters  r and  O' ' . 

(1.3.1)  f (T|6)  = — Tr_1  exp  (-T/6),  T-o 
r (r-1 ) ! 

Thus,  the  unconditional  probability  density  function  of  T is 

;co 

0 fr  (T|0)  <1  (01  X.  l)  dt> 


^ 1 > °- 

This  is  just  a scale  transform  of  the  inverted  beta  distribution  written 
in  this  form  to  show  its  similarity  with  (1.2.3). 

1.4  Existing  failure  data  (especially  field  data)  usually  do  not  exhibit 
any  of  the  two  features  discussed  above.  Usually  the  test  or  operational 
time  varies  from  equipment  to  equipment  of  the  same  type.  Thus,  the  data 
will  usually  be  of  the  form  (r^  T.),  i=l,  ....  n,  where  r^  is  the  number 
of  failures  of  the  ith  equipment  in  time  T^.  In  a test  situation,  it  is 
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feasible  to  control  either  r^  or  T^ , but  cost  considerations  recommend  the 
control  of  T. . Thus,  it  is  desirable  to  estimate  \ and  > in  this  more 
general  situation,  which  encompasses  the  situations  discussed  in  sections 
1.2  and  1.3  as  special  cases.  Schafer  et  al  [3]  present  a method  of 
estimation  akin  to  the  method  of  moments.  This  method  however  is  not 
applicable  if  a single  equipment  had  no  failures  at  all. 

1.5  In  this  report  we  present  a general  estimation  method  which  we  call 
The  General i zed  Maximum  Like! i hood  Method . A sufficient  condition  for  the 
existence  of  the  estimators  is  given.  In  the  case  of  fixed  time  data, 
it  is  shown  that  the  condition  is  also  necessary.  The  method  has  the 
advantage  of  being  usable  to  update  the  prior  when  new  data  become  available, 
e.g.,  from  reliability  demonstration  tests. 

If  the  data  used  for  the  estimation  of  the  prior  distribution  are  generated 
by  a planned  test,  the  estimability  condition  dictates  ways  of  choosing 
(control  1 ing)  ei ther  the  test  times  T..  or  the  number  of  failures  r.  in  such  a 
way  that  the  resulting  Generalized  Maximum  Likelihood  Equations  have  a 
solution,  i.e.,the  estimators  exist. 

In  the  case  of  fixed  time  data,  if  the  estimability  condition  is  violated, 
some  alternate  estimation  methods  are  presented. 


i 


i 

I 
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2.  THE  GENERALIZED  MAXIMUM  LIKELIHOOD  ESTIMATION  METHOD 
2.1  Ue  suppose  that  n identical  equipments  with  exponential  time-to- 
failure  distribution  are  tested  in  the  followinq  way:  the  ith  equipment 
is  tested  for  T ■ hours,  i = 1,  ....  n.  Let  r-  denote  the  number  of  failures 
of  the  ith  equipment.  We  assume  that  the  prior  distribution  of  the  MTTF  0 
is  given  by  (1.1.2).  Then,  the  unconditional  probability  function  of  r^  is 
given  by  (1.2.3),  i .e. , 


\ <',>  (i^r) 

The  Generalized  Likelihood  Function  of  the  sample  (r^,  T . ) , i 
is  defined  to  be 

l-,i,  (TKvO  (5r)" 


- 1 , . . . , n 


Just  as  in  the  classical  Maximum  Likelihood  Estimation  technique,  the  best 
explanation  of  the  data  (r^,  T.),  i = 1,  ...,  n is  provided  by  the  values 

A A 

(\,  y)  of  (\,  y)  at  which  the  function  L attains  its  maximum,  if  L has  a max- 
imum. As  usual,  in  order  to  maximize  L,  it  is  enough  to  maximize  its  natural 
logari thm. 

n / \+r.:-l  \ n v n T^ 

(2.1.3)  InL  * E In  ( ) + \ E In  — - — + Z r.  In  — 

i = l \ ri  / i = l T.+y  i = l T.+y 

In  order  to  obtain  the  critical  point  of  L,  we  have  to  solve  simultaneously 
the  General ized  Li kel ihood  Equations . 
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(2.1.4) 


3 InL  = o,  3 lnl  = o 


3X  3y 

which  in  our  case  become 


(2.1.5)  _3_1  nL  = Z / — + . . . + _J V " ln/'l+Ti)=0 

3X  i = l ^ X X+rrl  / i = l V T J 


(2.1.6)  _3_1  nL  = 

3Y  Y 


X 


i = l T-j+y  i = l 


Ti^ 


If  we  set  a.  = Z 1 the  above  equations  are  reduced  to 
J ri  2J 


(2.1.7)  Z 


J>1  X+j-1 


(2.1.8)  X = y Z 
i = l 


J,'"  (’  * T~) 


r . / n 

" / * 

Ti+y  / i = 


Ti 


1 Ti+y 


Since  X is  given  explicitly  in  terms  of  y by  (2.1.8),  we  can  substitute  in 
(2.1.7)  to  obtain  an  equation  in  y alone.  The  resulting  equation  can  be 

solved  numerically  (when  a solution  exists)  to  obtain  the  estimator  y and 

/\ 

then,  by  means  of  (2.1.8)  obtain  the  value  X . 

2.2  If  we  control  the  number  of  failures  r^  and  let  T.  be  random,  the 
distribution  of  Ti  is  given  by  (1.3.2).  It  is  immediate  that  the  new 
Generalized  Likelihood  Function  will  be  the  same  as  the  one  given  by 
(2.1 .2)  up  to  a factor 

n 

n r./T. 

1-1  1 1 

which  is  independent  of  the  parameters  X and  y.  Therefore,  the  resulting 
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Generalized  Likelihood  Equations  will  be  exactly  the  same  as  the  ones 
given  by  (2.1.7)  and  (2.1.8).  Thus,  the  Generalized  Maximum  Likelihood 
estimators  w i 11  h aye  t he_ _s_jni_e  form,  irrespective! y of  whether  we  control 
the  T^'s  or  the  r-'s  or  any  combination  of  them  (e.g.,  controlling  r^  for 
i = 1,  ....  k and  Tj  for  j = k + 1,  ....  n). 

3.  A SUFFICIENT  CONDITION  FOR  THE  EXISTENCE  OF  THE  GENERALIZED  MAXIMUM 
LIKELIHOOD  ESTIMATORS 

3.1  The  system  of  equations  (2.1.7)  and  (2.1.8)  does  not  always  have  a 

solution.  Although  we  could  produce  examples  of  actual  data  for  which  the 

Generalized  Maximum  Likelihood  estimators  do  not  exist,  for  simplicity's 

sake  we  resort  to  the  following  rather  contrived 

EXAMPLE  3.1.1  Let  n = 3,  r*  = o,  r = 1 , r.,  = 2,  T,  = T = T = T. 

i 2 J 12  3 

Then  the  equations  (2.1.7)  and  (2.1.8)  are  reduced  to 

(3.1.1)  2+1  = 31n  ( 1 + T \ 

X X+l  ^ Y J 

(3.1.2)  X = y/T 

whose  simultaneous  solution  calls  for  the  zero  of  the  function 


y(X)  = _i_+  - 31n  (1  + | ),  X > o. 


We  claim  that  actually  T (X)>o  for  all  X>o.  Indeed,  lim(X)=+°°  and  lim  H'(X)=o 

X-+o+  X-KO 

and  thus,  it  is  enough  to  show  that  V is  strictly  decreasing.  This  is 
so  since  the  derivative  of  T1  is  negative. 


r (X)  = - 


2 _L  + 3 

X2  (X+l)2  X( X+l ) 


X+2 

? (X+l)2 


< o 
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3.2  We  now  give  a sufficient  condition  for  the  solvability  of  the 
Generalized  Likelihood  Equ  ".ions  (2.1.7)  and  (2.1.8).  We  use  the  following 
notation: 


r = 


1 


I r.,  S 

1 r 


2 . 1 


i = l 


n 


_ l 2 

T = — 1 - Z T S 

n i=l  1 T 

Cov  (r,  T)  = 1 


J 2 

Z r - r 
n i=l  1 


r T2  - T , 


n i = 1 i 


n 


Z r.  T.-rT. 
n i=l  1 i 

We  shall  prove  the  following 

THEOREM  3.2.1  If 

(c)  2 r T cov  (r,  T)  < T (S  - r)  + r , 

r T 

then  the  Generalized  Likelihood  Equations  (2.1.7)  and  (2.1.8)  have  a solution. 
PROOF.  It  suffices  to  show  that  the  function  defined  by 


(3.2.1)  W(y)  = Z 


aj 


J>1  A(y)+j-l 


n t . 

Z 1 n ( 1 + — ) , y > o 

i=l  Y 


where 


X(y)  = Y 


n r • / n T • 

z -ll_  / Z — _L 

i=l  T^+y  j i"l  Tj+y 


has  a zero.  We  observe  that  a-|$o,  because,  otherwise,  all  aj  = o which 

implies  that  all  r^  = 0 and  the  condition  (C)  violated  (it  is  reduced  to 

o < o).  Since  X(y)  tends  to  zero  when  y tends  to  zero  and  since  a^>o  we 

get  lim  W (y)  = + 00  . 

Y -+o+ 
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Thus,  it  suffices  to  show  that  W(y)  is  negative  for  large  y.  To  this  end 


we  observe  that 


(3.2.2)  A(y)  -y[  Z r,  - z r^T,Y~^+  o(y'')]/["  t,  - " tV'*  ofY-')] 
i-i  i=l  / i=l  i=i  i 


— — y + r sT  - t 


+ o(l)  as  y + + 


Substituting  (3.2.2)  into  (3.2.1)  we  obtain 


(3.3.3) 


W(y)  = 


s «./(-£-  y + Ut  : T| 

J>1  J/  T 


cov  (r,  T)  + . 


j - 1 +o ( 1 ) ] - Z In  (1+  _1_) 
i=l  y 


_ 2 _ 

T -1  r r Sr  ■ T cov  (r,  T)  , 

v t ‘j  "t  = ~ ^ aj+  ~ — ^ (j"1)a1*]Y  + 

r & r T J>1  J r j>l  J 


-1  n _i  n 2 -2  -2 

* 0(.  )]  - Z Ty  + 1 E T.y  + o (y  ) 
i*!  2 i=l  1 


as  Y -»  ♦ 


2 

Since  » - nr  and  E (j-1)  a = 1 (S  +7-7) 

j>l  J j-1  J 2 r 


(3.3.3)  is  reduced  to 


W (y) 


n . _ _ 2 2 2 2 2 ? 

— ^ - 2 r T cov  (r,  T)  + T (S  - r")  + 7 S ) y"  + o(y"  ) 

r 


as  y + oo 


Because  of  the  condition  (C),  W(Y)  < o for  large  A and  the  theorem  isproved. 
3.3  We  observe  that  in  case  T1  = T for  all  1-1.  ....  n.  the  condition  (C) 


is  reduced  to 


(Ci)  7<S^ 
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which  is  exactly  the  condition  for  the  applicability  of  the  method  of 
moments  to  the  estimation  of  the  parameters  of  a negative  binomial  distri- 
bution. In  the  next  section  we  will  show  that  the  condition  (C^)  is  also 
necessary  for  Maximum  Likelihood  estimability  in  the  case  of  the  negative 
binomial  distribution. 

If,  on  the  other  hand,  rj  = r for  all  i = 1,  ....  n,  then  the  condition  (C) 


is  reduced  to 


T < fr 


which  is  exactly  the  condition  for  the  applicability  of  the  method  of 
moments  to  the  estimation  of  the  parameters  of  the  inverted  beta  distri- 
bution (1.3.2).  Of  course,  if  the  method  of  moments  is  applied  to  this 
situation,  the  value  of  \ will  always  be  greater  than  2,  because  this 
distribution  does  not  have  a second  moment  if  \ .<2. 

Another  remark  about  the  condition  (C)  is  in  order.  Since  the  values  of  Tj 
are  usually  large,  both  sides  of  the  inequality  become  rather  large.  It 
is,  therefore,  convenient  to  write  condition  (C)  in  the  following  equivalent 


form  (3.3.1) 


' -H,  (r-t  T’>: 


3.4  We  now  prove  the  following 

THEORLM  3.4.1  The  Maximum  Likelihood  estimators  of  the  parameters  of  the 
negative  binomial  distribution  (1.2.3)  exist  if  and  only  i f r < S 

r 

The  sufficient  part  of  this  theorem  is  contained  in  the  Theorem  3.2.1. 


For  the  necessity  part  we  need  the  following  lemma: 


g 


UMMA  3.4.2  For  dll  positive  Integers  n and  dll  positive  X we  have 
(3.4.1)  n 

>:  I , n 

J-l  (X+hV  \(*X+n-l ) 

The  Inequality  is  strict  if  n s 1. 

The  proof  of  the  lemma  is  inductive.  The  inequality  is  obviously  true  for 
n ■ l.  Assuming  it  true  for  n,  we  shall  prove  it  for  n + 1.  Indeed 
n+1  , n . 

'•  V . y .1  . . + I n n + 

j-1  (v+j-l)1’  j-l  (v+j-l)*’  (Hn)?  \(\>n-D 

J n+  If  n s n+  1 

\(\+n)  \(\+n-l)  (\+n)*  \ ( n ) 

This  completes  the  proof  of  the  lemma. 

PROOF  OF  Till  THEOREM:  We  need  only  prove  the  necessity  of  the  condition. 

Since  T ^ - T for  all  1 I n,  th  equations  (2.1.6)  and  (2.1.6)  are 

reduced  to 

n i « t 

(3  4.2)  y ( >...<■  ) - nln  ( 1 ♦ ) * o 

M \ Urr1  3 

(3.4.3)  T/y  - r/x. 

Obviously,  it  is  enough  to  prove  that  the  function  4 defined  by 
n , , 

3 ‘ ( \ ) y.  ( Or...*  ) - n 1 n ( 1 + r ) , \ s o 

1*1  X X+r.j-1  \ 

O 

has  no  zeros  if  r > S , If  all  r,  are  zero,  then  the  function  4 

r ' 

is  identically  zero,  the  log-likelihood  function  (2.1.3)  is  reduced  to 


U) 


nAl  n 


T+Y 


and  It  is  obvious  that  this  function  has  no  maximum  in  the  range  of  the 
parameters.  Thus,  we  may  assume  that  at  least  one  r j is  different  than 
zero.  Then 


lim  'F(A)  = + °°  and  lim  'f'(A)  = o. 

A -*•  o + A -*■  + °° 


Hence,  in  order  to  show  that  the  function  4*  has  no  zeros,  it  suffices  to 
show  chat  its  derivative  is  nonpositive.  By  using  the  inequality  (3.4.1) 
we  get 


(3.4.4)  •'  (A)  = - Z ( 


1 


+ ...  + 


1 


i =1 


) + 


nr 


(^+r1--l)^  A(A+r) 


< - 2 


nr 


i=l  A(A+rj-l)  A(A+r) 


We  now  use  the  convexity  of  the  function 

w (x)  = 1 / ( A- 1 + x),  x > 1 ; A > o 


with  the  rj's  in  the  numerator  of  the  summand  of  the  right  hand  side  of 
the  inequality  (3.4.4)  used  as  weights  and  the  r^ 's  in  the  denominator 
as  points  in  the  domain  of  w.  The  so-called  Jensen's  inequality  yields 

2 


ri 


EH 


n r 


i =1  A+rrl 


A-l  + Zr,  /zr1 


Ar  + S?  + r 2 - 
r 


Therefore,  going  back  to  (3.4.4)  we  get 

, - 2 

r (a)<  - _J d r 


nr 


A Ar+S^+r^-r  A(A+r) 
r 
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r 


n r (r  - S ) 
r 


\(X+r)  ( Xr  + S2  + r r ) 
r 


< o,  \ > o 


since 


2 2 

r 1 s . + r - f « l I (j-l)  a - o 

r r n j>1  J - 


This  completes  the  proof  of  the  theorem. 


4.  CONTROLLING  FOR  ESTIMABILITY 

4.1  If  we  let  n = 2m,  = T for  i =1,  ...  , m and  = kT  for  1 = m + 1, 

. .. , 2m,  k > 1 , then 


(4.1.1)  _ 

T = 


M t.  S* 


T , cov  (r,  T)  = 


k- 

2 


1 |"l  2m 

— T — — E 

„ m i =n 


By  substituting  (4.1.1)  into  (C)  we  obtain 


2m 


k+1 


(4.1.2)  2 r E r.  < 

m i=m  k-1 


r k+1 


a condition  independent  of  T.  This  kind  of  test  designing  enhances  the 
possibility  of  having  S2  > r and  the  condition  (C)  satisfied.  In  particular, 
if  k = 2,  the  condition  (4.1.2)  is  reduced  to 


(4.1.3)  2 r £ r.  < 3 (s^  -r)  + r 2 

m i=m  d 


4.2  In  the  more  expensive  case  of  controlling  the  number  of  failures  and 

_ 2 

letting  the  test  time  be  random,  we  can  always  assure  r < Sr  . We  consider 
the  following  design:  n = 2m,  ri  = r for  i =1,  ....  m and  ri  = kr  for 
1 « m + 1 , . . . , 2m,  k > 1 . Then 

2m 

a 1-1  I 1 

(4.2.1)  r = -ill—  r,  S‘ 


k-1 


S2  = ^_ili_^2r2,  cov  (r,  T)  = r 


“2m  “ 

L_  E Ti  -T 
m i=m  1 


and  the  condition  (C)  is  reduced  to 


(4.2.2)  2 T 

m 


T E T1  < T 2 fiJ-  r + 2 JilTI  + Jzl 
i=m  k+1  k-1  k+1 


*T  • 
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By  choosing  k = 2 and  r = 12,  or  k ■ 3 and  r = 7,  or  k = 4 and  r = 5,  or 
k = 5 and  r = 4,  this  condition  Is.  always  satisfied.  Of  course,  one  does  not 

have  to  go  to  such  extremes.  For  k = 2 and  r - 6 for  example,  the  condition 
will  usually  be  satisfied. 


5.  VARIABILITY  OF  THE  ESTIMATORS 


5.1  The  variance-covariance  matrix  of  the  Generalized  Maximum  Likelihood 


estimators  X and  y Is  given  by 


(5.1.1) 

“ 

var  (X) 

A A 

COV  (X,y) 

E 92  InL 

vF 

E 92  InL 

9X3y 

A A 

cov  (X,y) 

A 

var  (y) 

£ 92  InL 
9X9y 

E92  InL 
77 

9y  ; 

From  (2.1.5)  and  (2.1.6)  we  obtain 


(5.1.2)  ^lnl  ■ - A ^ + •••  * — .,) 


1=1 


(A«yl)2' 


(5.i.3)  ir.1nL  5 l A .It 


3X9y 

)2 
3y 


y 1=1  T-PT 


lfj"L  = - _X_  E Ti^+  V + " 
(5.1.4)  .X  ^ LI  -^JT" 


5.2  Assuming  the  T^  non-random  and  the  r^  random. 


(X+j-1) 


2 


(T,  +Y)2’ 
we  obtain 


(5.2.1)  E _J2  InL  = - E F ( Tj  ; \) 
~9X2  1s1  Ti  + Y 
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where 


X+j-1 


(5.2.2)  F ( p ; X ) = I pj/j2  ( j ) 

j-1 


<5-2-3'  E »2  1„L  - 1 Z T. 


3A9y 


y 1=1 


Ti 


/ c o /IN  r ^ 1 nL  = - X Z 

(5.2.4)  E— r — 2 i=1 

3y  y 


Ti 


Ti 


By  substituting  (5.2.1),  (5.2.3)  and  (5.2.4)  in  (5.1.1)  we  get 


(5.2.5)  var  (A)  = A/  [A  Z F (Jj ;A)  - IT.  ]; 

1=1  T.  +y  T~ +y 

a n n y Ti  ^ Ti  ^ Ti 

(5.2.6)  var  (y)  = y2  I F ( Ti  *A)/  E — [A  I F (— L-  ;A)  - I -1—  ] 

i =1  T +Y  jr  i=1  Ti+T'  i=1  Ti+Y  i=l  Ti+V 

(5.2.7)  cov  (A ,y)  = y/[A  I F(]jt ;A)  - I Ti  ]. 

i-1  T.  +y  i =1  +y 


5.3  By  assuming  T.  random  and  r.  non-random,  we  obtain 


(5.3.1)  E 3 


lnL= . I aj  o . 
3F"  j>1  ( A+  j-1r  ’ 


(5.3.2)  E 8 InL  = - I ri 


8A9y 


1=1 


A +r. 
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(5.3.3) 


E 3 


1 nL 


3Y 


1-1 


X + r1  + 1 


Substituting  into  (5.1.1),  we  get 


(5.3.4)  var  (X)  = X E i 


1=1 


X + r . , +1 

a . 


/ A; 


(5.3.5)  var  (>)  = / >,(^2  / A; 


(5.3.6)  cov  (X.y)  = Y £ ri  / A , 

1=1  X+  r.. 


where 


(5.3.7)  A=  X E “j 


r . 


H O 

- (I  _J )2. 

i=l  X+r^ 


The  case  of  mixed  controls  can  be  handled  similarly.  In  order  to  estimate 
these  variances  and  covariances,  we  substitute  in  the  above  formulas  the 

A A 

estimated  values  X and  y instead  of  X and  y 


6.  NUMERICAL  CONSIDERATIONS 

6.1  In  writing  up  a computer  program  for  the  numerical  solution  of  the 
Generalized  Maximum  Likelihood  Equations  (2.1.7)  and  (2.1.8)  the  following 
observations  should  be  taken  into  account.  The  solution  of  these  equations 
is  reduced  to  finding  the  zero  of  the  function  W defined  by  (3.2.1).  The 
shape  of  this  function  is  given  by  figure  1. 
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If  the  Newton-Raphson  Iterative  method  is  employed,  care  should  be  exercised 
In  choosing  the  initial  value.  If  the  initial  value  is  greater  than  the 
minimum  m of  the  function  W,  the  Newton-Raphson  process  will  diverge  to 
infinity  (initial  value  Y-j  in  Figure  1),  whereas,  if  the  initial  value  is 
less  than  m but  near  m,  the  first  iteration  will  produce  a negative  value 
for  y (initial  value  yQ  in  Figure  1).  Therefore,  an  initial  value,  at  which 
W is  positive,  should  be  chosen,  if  the  Newton-Raphson  method  is  to  be  used. 
Because  of  the  complexity  of  the  derivation  of  W and  since  only  nearest 
integer  accuracy  is  required  for  y,  some  slower  converging  interpol ati ve 
method  may  be  more  suitable. 


7.  ALTERNATE  ESTIMATORS  (NEGATIVE  BINOMIAL) 

7.1  In  the  case  of  fixed  time  testing  (l.e.  T.  = T for  all  i - 1,  ....  n), 
we  saw  that  the  distribution  of  the  number  of  failures  is  the  negative 
binomial  given  by  (1.2.3)  and  that  the  Maximum  Likelihood  estimators  exist 
1 f and  only  i f 

(7.1.1)  7<S? 

r 

by  the  theorem  3.4.1.  Since  the  solution  of  the  Likelihood  Equations 
require  numerical  techniques,  several  alternate  methods  are  often  used. 

Among  them,  the  method  of  moments  is  the  most  popular.  It  yields 

(7.1.2)  A ? A , 

' = r / (S^-r)  . > = r T / (S^-r) 

and  it  is  highly  efficient  in  a wide  range  of  the  parameters.  Of  course, 
this  method  is  usable  if  and  only  if  (7.1.1)  is  satisfied. 

Other  simple  methods  are: 

(A)  Matching  first  moments  and  first  frequencies  (the  zero  class  of  the 
sample  with  the  expected  number  in  the  zero  class).  The  resulting  equations 
are 


(7.1.3)  \T  = r>,  n^  = (V/(r  +x))X  > 

n 

where  nQ  is  the  zero  class  of  the  sample.  It  is  not  hard  to  prove  that  this 
method  is  usable  if  and  only  if 

(7.1 .4)  r>  In  (n/no). 
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(B)  Matching  first  moments  and  the  ratio  of  the  first  two  frequencies. 

The  resulting  estimators  are 

A A 

(7.1.5)  A = nQr/(n1r-no)  , y=  nQT/(n1r-no) . 

Obviously  this  method  is  usable  if  and  only  if 

(7.1.6)  r > n0/n^  . 

The  efficiency  of  these  two  methods  has  been  investigated  by  Katti  and 
Gurland  [5].  They  found  that  there  are  ranges  of  the  parameters,  where 
these  methods  are  superior  to  the  method  of  moments. 

7.2  Another  estimation  investigated  by  Katti  and  Gurland  [5]  is  the  minimum 
chi-squared  method.  It  is  a highly  efficient  and  rather  complicated  method 
for  which  numerical  techniques  are  required.  We  did  not  attempt  to  find 
necessary  or  sufficient  conditions  for  the  applicability  of  this  method. 

8.  CONCLUSIONS  AND  RECOMMENDATIONS 

In  this  report  we  presented  a method  for  estimating  the  shape  and  scale 
parameters  of  an  inverted  gamma  prior  distribution  of  the  mean-time-to- 
failure  for  equipment  having  exponential  time-to-fail ure  distribution.  All 
sorts  of  existing  failure  data  on  the  equipment  in  question  are  usable 
provided  a certain  sufficient  condition  is  satisfied.  Further,  the  method 
can  be  used  to  update  the  prior  when  new  failure  data  become  available. 

This  periodic  updating  will  give  rise  to  a solid  prior  which  can  confidently 
be  used  in  Reliability  Demonstration. 

It  is  recommended  that  a computer  program  be  written  to  solve  the 
Generalized  Likelihood  Equations  that  define  these  estimators  and  to  compute 


their  variance-covariance  matrix.  To  this  end.  the  recommendations  pot 

I.  Section  « shooid  be  taken  into  accoont.  The  pro, ram  can  then  be 
used  for  the  periodic  updating  of  the  prior  distribution. 

Examples  of  the  application  of  the  theor,  developed  in  this  report  are 
presented  in  the  Appendix. 
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APPENDIX 

In  this  appendix  we  present  some  examples  of  failure  data,  for  which  the 
theory  developed  in  this  report  is  applied. 

EXAMPLE  1 . Three  pieces  of  the  equipment  AN/ASK-6  were  tested  with  the 
following  results  (time  measured  in  hours): 

r]  = 4,  T]  = 1961;  r2  = 2,  T?  = 1814;  r3  = 1 , T3  = 1890 

thus  r = 7/3,  T = 5665/3  and  the  right  hand  side  of  (3.3.1)  is  1.44271, 
whereas  the  left  hand  side  is  2.3333.  Therefore,  the  condition  (C)  is  noj 
satisfied. 


EXAMPLE  2.  Equipment:  AN/ASN-108,  n=3,  r]  = 3,  T]  = 1522;  r?  * 0. 
T,  = 1725;  r3  = 1,  T3  = 997 

It  is  easily  verified  that  the  condition  (3.3.1)  is  satisfied 

(1.333  < 1.699).  Obviously  a-j  = 2,  a2  = a3  = 1 and  = o for  -4. 

and  the  generalized  maximum  likelihood  equations  (2.1.7)  and  (2.1  8 ■'** 

(A1  ) W(y)  = 2 + 1 + 1 - 1 n ( 1 + 1522)  - 1 n(  1 ♦ 1 725,  - 

A(y)  A(y)+1  a(y)+2  y 

1 n(l  + 997)  = o 

Y 


where 
(A2 ) 


A(y) 


( 


3 + 1 

1 522  + y 957 


1 \/Yl522 

997  + yjj  ^1  522  + 


+1725  + 997 

7 v ■,  wtt 


*0 


In  order  to  find  the  estimators  A and  y , we  start  with  a moderate 
initial  value  yQ  = 1000  and  find  W(1000)  = .267688.  Since  W(yq)  so 
we  try  y-j  = 2000  and  find  W(2000)  = .039124.  By  linear  extrapolation  we 
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OWN 


try  the  next  value: 


>2  - >1 


W(>,  ) - 2,171 .17 


'1  “'o 

wTvjT^  w7V0T 

and  find  W(\0)  * .028774.  Linear  extrapolation  produces 


1 


\ 


1 y ~ 1 1 

Y,  « 1, W(yJ  - 2.647.04 

w(>2)  - w(> 1 ) 

with  W(>  ,)  = .011590.  Continuing  this  way,  we  find  ^ 2968 

W(>4)  = .005417,  >5  = 3250,  W(>5)  * .001931,  - 3406,  W(>fi)  -.000525 

- 3464,  W(>  7)  = 76. 1807. 10'6,  - 3474  , W(>  ) « 24.1  77.10"7 

W ( 34 75)  ^0  . Actually  > = 3474  \ = 3.332. 

In  this  example,  we  approximated  > by  staying  on  the  left  side  of  its 
true  value  and  by  using  linear  extrapolation.  Of  course  we  could  have  started 
with,  say  >o  = 2000,  = 4000,  observe  that  W(yq)  No  and  W(\j)  - and  use 

linear  interpolation  to  find  > ^ and  continue  this  way.  Of  course  in  this 
example,  the  number  of  the  equipments  is  too  small  to  rely  on  the  estimators. 
It  is  only  presented  for  the  purpose  of  demonstrating  the  numerical  technique. 

A \ 

The  variability  of  \ and  ^ can  be  computed  from  the  formulas  (5.2.5),  (5.2.6) 
and  (5.2.7)  if  we  assume  that  the  T.'s  are  controlled  and  the  r.'s  are  random. 
This  assumption  was  not  exactly  met  in  this  experiment,  but  it  is  more 
realistic  than  the  other  way  around.  We  find 

A 

standard  deviation  of  \ - 10.27 
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standard  deviation  of  y 

/V  A 

covariance  of  (A,y) 
correlation  coefficient 


= 10,903.56 
= 109,985.92 
of  (A.y)  = .9821 


Thus,  the  variability  of  the  estimators  is  extremely  high,  which  is  expected 
in  view  of  the  low  value  of  the  number  of  equipments  tested. 


EXAMPLE  3.  The  following  field  failure  data  of  a Ground  Electronic 
Oigital  Processing  System  were  collected: 


r 

T 

r 

T 

1 . 

3 

5068.6 

17. 

6 

6435.2 

2. 

15 

7486.2 

18. 

7 

4624.5 

3. 

10 

7587.4 

19. 

14 

5327.0 

4. 

5 

4978.4 

20. 

11 

7486.2 

5. 

7 

6000.0 

21  . 

3 

6271  .7 

6. 

5 

5187.5 

22. 

9 

6934.5 

7. 

3 

7808.5 

23. 

16 

7114.8 

8. 

3 

2246.4 

24. 

9 

7626.1 

9. 

7 

4735.2 

25. 

7 

4372.2 

10. 

13 

7670.9 

26. 

10 

5409.2 

11  . 

3 

4320.2 

27. 

8 

5617.6 

12. 

3 

3599.1 

28. 

4 

2844.8 

13. 

6 

7865.8 

29. 

10 

1976.3 

14. 

1 

1941 .5 

30. 

7 

1987.3 

15. 

21 

7273.6 

31  . 

3 

4952.5 

16. 

17 

6891  .8 

The 

condition 

(3.3.1) 

is  satisfied. 

The 

values  of 

the  a's  an 

a2 

= 30 , = 

30, 

= 23,  ctg  = 22,  oig  = 

20,  a-j  = 18,  dg  = 13 

a10 

= 10,  a1 1 

S 2 “ 6,  ^1  2 ~ ® » 

al  4 

= 5*  “15  = 

4,  alfi  - : 

al  8 

" a19  = “ 

20  ~ “21 

= 1. 

By 

starting  with  the 

values  ^ = 3000 

and 

Y-|  = 5000, 

we  obtain 

v = 3000  \ 

= 4.4313261 

= 

.41345261 

0 

0 

0 

Y-j  = 5000 

= 7.3422691 

W1  = 

> -.07330695 
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v2  - Yi  - W1  (Yl  - Y0)  /(w,  - W0)  « 4700 

X2  = 6.9017913  W2  = -.0328700 
Y3  = Y2  - W2  (y2  - Yj) / (W2  - W,)  * 4456 

X3  = 6.5514034  W3  = -.0128367 
Y4  = Y3  - W3  (y3  - Y 2 ) / (W3  - W2)  - 4300 
Y4  « 6.3297434  W4  = -.0040735 

Since  calculations  on  a pocket  calculator  are  very  tedious,  we  stop  by  just 
extrapolating  once  more.  Thus, 

Y - Y4  - W4  (y4  - Y3) ] (W4  - W3)  = 4222 

X = X4  - W4  (X4  - X3) J (W4  - W3)  - 6.227 

The  variability  of  these  estimators  is  relatively  low.  For  example,  the 
s.d.  of  X computed  by  means  of  the  formula  (5.2.5)  is  approximately  2.918. 
Large  test  times  and  large  number  of  equipments  tested  naturally  decrease 
the  variance  of  the  estimators. 
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MISSION 

of 

Rome  Air  Development  Center 

RAVC  plant*  and  execu. tes  research,  development,  test  and 
selected  acquisition  programs  In  support  of  Command,  Control 
Communications  and  Intelligence  IC3I)  activities.  Technical 
and  engineering  support  within  areas  of  technical  competence 
Is  provided  to  ESV  Program  Offices  IPOs)  and  other  ESV 
elements.  The  principal  technical  mission  areas  are 
communications , electromagnetic  guidance  and  control,  sur- 
veillance of  ground  and  aerospace  objects,  Intelligence  data 
collection  and  handling,  Information  system  technology, 
Ionospheric  propagation,  solid  state  sciences,  microwave 
physics  and  electronic  reliability,  maintainability  and 
compatibility. 


